Non— universality of elastic exponents in random bond— bending networks 
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We numerically investigate the rigidity percolation transition in two-dimensional flexible, random 
rod networks with freely rotating cross-links. Near the transition, networks are dominated by 
bending modes and the elastic modulii vanish with an exponent exponent / = 3.0 ± 0.2, in contrast 
with central force percolation which shares the same geometric exponents. This indicates that 
universality for geometric quantities does not imply universality for elastic ones. The implications 
of this result for actin-fiber networks is discussed. 
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In contrast with ordered crystals, disordered materi- 
als typically exhibit a complex connection between their 
collective mechanical properties and the underlying in- 
teractions between their constituent elements . In this 
Rapid Communication, we study this connection between 
effective elastic theories of a flexible rod network in two 
dimensions and the fundamental elastic properties of in- 
dividual rods, along with statistical measures of the net- 
work connectivity. We focus on cross-link densities to 
explore the break-down of collective rigidity of the net- 
work. This vanishing of the static shear modulus in re- 
lated systems 0, H, S S IE 13 has been termed the rigid- 
ity percolation transition in analogy to the better under- 
stood scalar, or conductivity, percolation transition |lflj . 
Previously we investigated the internal deformation 
field of a 2d flexible rod network away from the criti- 
cal point of rigidity percolation, and found a cross-over 
from affine to nonaffine deformation upon changing the 
mean density between cross links in the network and/or 
the inherent flexibility of individual rods, which wasl also 
independently corroborated by another group 3j . 

The rigidity percolation transition has been previously 
studied numerically in a number of model, disordered 
systems with variations in both the nature of the links 
(bonds) and the disorder of the network. While the rigid- 
ity transition is first order on the Bethe lattice and in 
random networks of infinitely rigid rods , in all other 
random network models it is continuous and character- 
ized by a diverging length scale over which the material 
acquires a finite, static shear modulus. From the out- 
set, it is important to distinguish between two different 
classes of physical observables when discussing the long 
length scale physics at the transition: (i) geometric quan- 
tities which describe the fractal structure of the percolat- 
ing rigid cluster at the transition; and (ii) elastic prop- 
erties of the material near the percolation point. In the 
first category are the exponents v and /3, which describe 
the divergence of the correlation length as one approaches 
the transition and the probability that a link will be part 
of the percolating cluster, respectively. In the second 
category is /, which describes the power law dependence 
of the system's elastic moduli upon the approach to the 
percolation transition from the rigid phase (see Eq. ^ . 

In lattice models, the control parameter, as for scalar 



percolation, is the probability of the presence of a link in 
the diluted lattice. It appears that the nature of disor- 
der is a relevant variable, since for purely central forces 
between network nodes, site and bond disorder exhibit 
quantitatively different scaling regimes near the critical 
point 0, In addition, a third universality class has 
been postulated in the so-called "bond-bending" model, 
where bending lattice edges and rotation at vertices cost 
energy ^31 • This suggests that the introduction of a 
bending modulus in the model is a relevant perturbation 
at the transition. 

Pertinent to the current work are the simulations of 
Latva-Kokko et al. These investigations, depart- 

ing from previous work, turned to off-lattice simulations. 
They constructed two-dimensional random rod networks 
and applied a topological approach to investigate the 
rigidity percolation point in random rod networks with 
Hookean central forces and a bending modulus. The in- 
troduction of a bending modulus is vital at this stage 
since the random spring network that is not prestressed 
has zero-frequency deformation modes at any finite den- 
sity of links, and is thus always nonrigid (14] . Latva- 
Kokko et al. considered two variants of their model of 
fiexible rods distinguished by the constraint forces im- 
posed at a cross link: one with cross links that fix the an- 
gle between intersecting rods by applying local constraint 
torques, and freely rotating bonds at cross links in the 
other. In both cases they found that the geometric expo- 
nents are consistent with those of the rigidity percolation 
transition in a diluted lattice with central forces. This re- 
sult suggests that the introduction of bending forces are 
not relevant perturbations at the critical point in the fol- 
lowing restricted sense: the exponents associated with 
length scale and geometry of the spanning, rigid cluster 
appear to be universal. We use this result later. Their 
approach, however did not allow them to investigate the 
scaling exponent for the elastic quantities of their model. 

This Rapid Communication presents a study of the 
elastic properties of the rod network identical to the 
model system studied by Latva-Kokko et al.. We, how- 
ever, concentrate on the mechanical aspects of the net- 
work near the transition. It is important to note that 
this work relates to models of semi-flexible polymer net- 
works with freely rotating bonds Tsl. This detail of the 
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model distinguishes it from its other close antecedent, the 
bond-bending model of S. Feng et al. That lattice-based 
model has been particularly well described by both nu- 
merical simulation and real-space renormalization group 
techniques j^]. Our model differs in that the two fila- 
ments that cross at each node contribute independently 
to that node's bending energy; there is no energy cost 
for relative rotation between rods. Previous lattice cal- 
culations suggest that the scaling of the elastic constants 
near the transition depend on such details of the network; 
with this in mind, we seek herein to explore the rigidity 
collapse of sparse actin networks. 

The principal result of this communication is that the 
mechanical properties of the flexible rod network at the 
rigidity percolation transition are distinct from previ- 
ously investigated models of either lattice-based bond- 
bending networks or central force networks. This is true 
despite the fact that the scaling of the size of the per- 
colating cluster and its fractal geometry (as determined 
by Latva-Kokko et al. 0]) are identical to that of central 
force networks. This point highlights the physical in- 
dependence of geometric/topological exponents and the 
elastic exponents of the network. We propose that whilst 
some degree of universality in the elastic properties of dis- 
ordered systems may exist, such universality classes are 
smaller and more numerous than for geometric proper- 
ties; that is, models with the same geometric exponents 
V and [3 may have distinct elastic exponents /. Our 
model system below provides a concrete example of this 
claim. We suspect that the underlying cause of this non- 
universality is that while the rigidity transition itself de- 
pends on large scale topological properties of the network 
that are insensitive to the local details of stress transmis- 
sion, the moduli of the rigid network near the transition 
are sensitively determined by only a vanishingly small set 
of paths for stress propagation. The moduli thus depend 
on local details of how stress is transmitted through the 
small set network nodes along the weak points of this 
path. Different force laws at e.g. cross links can there- 
fore have profound effects on the way in which the moduli 
vanish at the transition. 

We numerically evaluated the elastic moduli of random 
networks of rods at T = 0, with an energy 871 per unit 
length 5s given by 

11 ~ 2\J's) ^ 2\J's) ^' 

The first term on the right hand side describes an elastic 
restoring force for changes in relative length of the rods, 
61 /Ss, with a spring modulus fi. Bending by an angle 69 
incurs an energy cost given by the second term, where the 
bending modulus k is the same as in the wormlike chain 
model (1^]. Cross linked rods are coupled by imposing 
the same coordinates at intersections. 

The simulation method is described in detail else- 
where [T^ (where we also discuss behavior away from 
the transition), but in brief: rods of length L are de- 
posited with random orientation and position into a two- 



dimensional shear cell of dimensions W x W. Each in- 
tersection is identified as a cross link, the mean distance 
between which (as measured along a rod) is denoted Ic , 
so that the mean number of crosslinks per rod is L/lc — 1. 
Deposition continues until the required cross linking den- 
sity L/lc has been reached. The system Hamiltonian 
is constructed from (Q, and the mechanical equilibrium 
configuration found by the preconditioned conjugate gra- 
dient method under the constraint of an applied shear or 
uniaxial strain. An example is given in Fig. ^ 




FIG. 1: ( Color online) The energy density for a network at the 
transition L/lc ~ 5.933 under an applied shear strain. The 
line thickness is proportional to the logarithm of the energy 
density per unit length. Apparent stressed 'dangling' ends are 
artefacts of numerical noise and make no significant difference 
to the measured quantities discussed below. The calibration 
bar shows what proportion of the energy is due to stretching. 
For clarity, a small W = 7^L shear cell is shown. 

There are four lengths in the problem: the system 
size and rod length W, L respectively; the mean dis- 
tance between cross links Ic and a length characteriz- 
ing the flexibility of the rods, Ih = \fnf\i. We take 
L/^c as our dimensionless measure of density rather than 
q = NL'^ with N the number of rods per unit area, as 
in 0, It is straightforward to convert between 

the two measures using the expression derived in the 
Appendix. Rigidity percolation occurs at (L//c) trans ~ 
5.933 6J. For Ljl^ slightly above this critical point, where 
£ — {L/lc)/{L/lc)trans — 1 0+, we writc the shear mod- 
ulus G as 

with a similar expression for the Young's Modulus, Y . 
For e > 0, and sufficiently large system sizes W , both G 
and Y depend only on e and l^/L. Fig.[21shows G and Y 



TABLE I: Quoted / for 2d uncorrelated networks 



Model 


/ 


Ref. 


Site percolation 


1.35 ±0.06 




Central force percolation 


3.57 ±0.3 




Bond-bending model 


3.96 ±0.04 


[i^, 20] 


This work 


3.0 ±0.2 


n/a 



versus e for systems where this large-Ty limit has been 
reached. Both moduli vanish continuously as £ — > 0"*", 
confirming both the nature and location of the transition 
as claimed in to within our data resolution [l8j |. 

Fig. shows that GL^/k is independent of l^/L for 
sufficiently small e, i.e. G is independent of n near the 
transition. This suggests that the transition is dominated 
by bending modes, and we can infer that G would vanish 
if there were only central force terms in agreement with 
[l4| . We can confirm the dominance of bending modes at 
the transition by measuring the proportion of stretching 
energy for increasing W. As shown in inset of this figure, 
this fraction is always small, <4%, and may vanish as 
— > oo, as suggested by the figure and in agreement 
with both 01 and our observation that G oc k. 

By fitting the data to the functional form 

G^Aef'{l + Be) (3) 

we find / = 3.0 ± 0.2, consistent with 3, where the 
error estimate gives the spread of values when the fit- 
ting procedure is repeated over different data subsets. 
B is always 0{1), suggesting is a 'sensible' choice. 
This is compared to known values in Tabled The only 
possible equality is with central force bond percolation; 
however, since the modulus near the transition is domi- 
nated by bending forces, we do not believe that this near 
agreement is meaningful. The fit for G can be made to 
agree with the Young's modulus data by simply rescal- 
ing A ^ A' with / and B fixed, as shown in Fig. |21 for 
Ih/L = 0.006. A' /A 2.7 ±0.2, so that the Poisson ratio 
a = Y/2G — 1 = 0.35 ±0.1 at the transition, consistent 
with the claimed 'universal' value of i ,19,] . 

Independent confirmation of / can be found by placing 
the system at the criticalpoint e = and considering the 
variation of G with |3j B • Scaling arguments sug- 
gest that G ~ l'^ for sufficiently large W ^ where v is 
the exponent describing the divergence of the correlation 
length, which is known for these networks to be in the 
range of = 1.17 ± 0.02 We attempted to fit our 
data to the same functional form used in ^9] , but the pre- 
cision with which we can extract the exponent using this 
approach is not sufficient to distinguish our results from 
either central force or bond-bending models. The num- 
ber of networks for each attempted W ranged from 400 
to 2000, comparable to lattice models which give much 
smaller errors. Our lack of precision is likely due to the 
random nature of our networks, in which the distance be- 
tween crosslinks varies continuously down to zero so that 
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FIG. 2: Dimensionless shear modulus (lower data set) and 
Young's modulus (upper set) for l\ylL — 0.006 versus e on 
log-log axes. The best-fit of to the data for G is plotted, 
and then shifted vertically to show agreement with Y. 



the strength of coupling between connected nodes can 
vary greatly, introducing a significant additional noise 
source into the problem. 
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FIG. 3: The dimensionless shear modulus vs. e for two differ- 
ent Ih/L on log-log axes. For clarity the error bars are not 
plotted in this figure (see Fig. I^J, but are no larger than the 
symbols. The solid line is a best fit to (|3Jl for l\^/L — 0.006 
and £ < 1, with each point weighted by its fractional devi- 
ation from the fit. (Inset) The fraction of stretching energy 
to the total, p, as a function of system size for Ih/L = 0.006. 
The line is a best-fit to the scaling form in 0]. 

In summary, we have shown that the elastic moduli of 
fiexible random networks with freely rotating crosslinks 
scale as near to rigidity percolation. The modulus it- 
self is dominated by bending modes near the transition 
and vanishes with an exponent / = 3.0 ±0.2. This expo- 
nent differs from previous reported exponents for models 
that include bond bending. Our introduction of bend- 
ing forces, motivated as it was by modelling cross linked 
F-actin networks, differs in detail from previous work. 
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Our results suggest that, in contrast to the geometric ex- 
ponents previous reported for this system, the modulus 
exponent is more highly model dependent. To address 
the experimental implications of these results for physi- 
cal actin networks, we note that the continuous transi- 
tion with true scale invariant behavior exists at only zero 
temperature. However, zero temperature phase transi- 
tions can have experimental consequences at finite tem- 
perature as long as it is possible to experimentally 
access the critical region. We expect that in sparse actin 
networks the algebraic decay of the static shear modulus 
will be cut off only by the entropic elasticity of the net- 
work |22|. a small quantity for a stiff polymer network. 
The experimental observation of this predicted decay, a 
signature of the T = transition at finite temperature, 
may be possible, but such quantitative predictions await 
a more complete description of the critical regime. 
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Appendix: Here we derive the relationship between 
L/lc and the quantity q — NL'^, where N is the number 
of rods per unit area, as used in 0, |E S Q| • Suppose a 
rod of length L lies along the x-axis oi a W x W box, 
where W ^ L so that the box shape should not matter. 
The probability of a second rod, deposited with orienta- 
tion 9 to the X-axis and random position, intersecting the 
X-axis is (L/W) sin 9. Hence the probability of intersect- 
ing the first rod is {L/W)^ sin9, which, when uniformly 



averaged over all 9 G [0, tt], becomes p = 2L'^ / {ttW"^). 

Since the rods are deposited at random, p is indepen- 
dent of how many other crosslinks each rod has. The 
mean number of intersections per rod is therefore simply 
p times the total number of rods in the system NW'^, 
i.e. pNW^ = 2L^N/'K = 2q/'K (assuming NW"^ > 1). 
The distribution P„ of the number n of cross links for 
any given rod is Poisson with a mean 2q/Tr. If each in- 
tersection is imagined to 'cut' the rod into n+1 line seg- 
ments, then the mean length of all such segments is just 
L/{{n) -I- 1) = L/{2q/Ti -f 1). However, Ic is defined be- 
tween crosslinks, not between crosslinks and the ends of 
the rods (which are dangling and thus do not contribute 
to the stress). In fact, only n — 1 segments contribute 
to Ic . Thus to evaluate Ic , we must average the quantity 
L/ {n + 1) over all valid P„ with a weighting ti — 1, and 
normalize accordingly, i.e. 



(4) 



Using the shorthand A = Zq/tt and the standard proper- 
ties of the Poisson distribution, we find: 



A + e~ 



1 



l+C-^ + f(c->-l)' 

which is a monotonic increasing function of A. Note 
that L/lc ^ 2(7/7r as g — > oo, and also that L/lc ~> 3 
as g — > 0+, as it should since (in this limit) the domi- 
nant contributions to Ic will come from rods with n = 2 
crosslinks, for which Ic is indeed L/3. 
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